Computerized Detection of Breast Cancer on Digital Tomosynthesis Mamograms

ABSTRACT

A method for using computer-aided diagnosis (CAD) for digital tomosynthesis mammograms (DTM) including retrieving a DTM image file having a plurality of DTM image slices; applying a three-dimensional analysis to the DTM image file to detect lesion candidates; identifying a volume of interest and locating its center; segmenting the volume of interest by a three dimensional method; extracting one or more object characteristics from the object corresponding to the volume of interest; and determining if the object corresponding to the volume of interest is a breast lesion or normal breast tissue.

CROSS-REFERENCE TO RELATED PATENT APPLICATIONS

This application is a continuation of U.S. patent application Ser. No.11/350,301, which was filed on Feb. 1, 2006 and claims the benefit under35 U.S.C. §119(e) of U.S. Provisional Application Ser. No. 60/650,923filed Feb. 8, 2005, the disclosure of which is hereby incorporated byreference in its entirety for all purposes.

FIELD OF TECHNOLOGY

This relates generally to breast cancer detection and, moreparticularly, to a system and method for using computer-aided diagnosis(CAD) for Digital Tomosynthesis Mammograms (DTMs) for detection andcharacterization of breast lesions.

DESCRIPTION OF THE RELATED ART

Cancer is a serious and pervasive medical condition that has garneredmuch attention in the past 50 years. As a result there has and continuesto be significant effort in the medical and scientific communities toreduce deaths resulting from cancer. Mammography is the mostcost-effective screening method for early detection of breast cancer.However, mammographic sensitivity is often limited by the presence ofoverlapping dense fibroglandular tissue in the breast. The denseparenchyma reduces the conspicuity of the abnormalities, whichconstitutes one of the main causes of missed breast cancer. The adventof full field digital detectors offers opportunities to develop newtechniques for improved imaging of dense breasts such as digitaltomosynthesis, stereomammography, and breast computed tomography. Thesenew techniques are still under development and their potential impactson breast cancer detection remain to be investigated.

Digital tomosynthesis is based on the same principle as conventionaltomography in diagnostic radiology which uses a screen-film system asthe image receptor for imaging body parts at selected depths. Inconventional tomography, a series of projection exposure is accumulatedon the same film when the x-ray source is moved about a fulcrum whilethe screen-film system is moved in the opposite direction. A drawback ofconventional tomography is that each tomogram can only image one planeat a selected depth in relatively sharp focus. If the exact depth ofinterest is not known in advance or the abnormality encompasses a rangeof depths, a tomogram at each depth will have to be acquired at separateimaging, thus costing additional dose and examination time.

With a digital detector, the series of projection exposure is readout asseparate projection views at the different x-ray source locations.Tomographic slices focused at any depths of the imaged volume can thenbe generated with digital reconstruction techniques from the same seriesof projection images. Because of the wide dynamic range and the linearresponse of the digital detector, each of the projection images can beacquired with a fraction of the x-ray exposure used for a regularprojection radiograph. The total dose required for digital tomosynthesisimaging may be kept at nearly the same or only slightly higher than thatof a regular radiograph. Properly designed digital reconstructiontechniques provide the additional advantage that the depth resolution oftomosynthesis is generally much higher than that of conventionaltomography. Digital tomosynthesis thus makes tomography more practicalto be applied to breast imaging in terms of radiation dose, examinationtime, and spatial resolution.

Digital breast tomosynthesis mammography is one of the promising methodsthat may reduce the camouflaging effects of dense tissue and improvemammographic sensitivity for breast cancer detection in dense breasts.DTM may also improve the accuracy in differentiation of malignant andbenign breast lesions because the reduction of overlapping tissue allowsthe features of the lesions to be visualized or analyzed more reliably.Several research groups are developing digital tomosynthesis methods forreconstruction of tomographic slices from the series of projectionimages. A study is underway to compare DTM with conventional mammogramsin breast cancer detection.

SUMMARY

Digital tomosynthesis mammography (DTM) is a promising new modality thathas the potential to improve breast cancer detection, especially indense breasts. However, the number of slices per breast may range from30 to over 80, thus increasing the time required for interpretation. Acomputer-aided diagnosis (CAD) method and a system is disclosed thatprovides a second opinion to radiologists during their interpretation ofDTMs. The disclosed CAD system includes retrieving a DTM image filehaving a plurality of DIM image slices; applying a three-dimensionalanalysis to the DTM image file to detect lesion candidates; identifyinga volume of interest and locating its center; segmenting the volume ofinterest by a three dimensional method; extracting one or more objectcharacteristics from the object corresponding to the volume of interest;and determining if the object corresponding to the volume of interest isa breast lesion or normal breast tissue.

BRIEF DESCRIPTION OF DRAWINGS

FIG. 1 is a block diagram of a computer aided diagnostic system that canbe used to perform breast cancer screening and diagnosis based on aseries of DTM images using one or more exams from a given patient;

FIG. 2 is a graph illustrating a distribution of a longest diameter of23 masses and 3 areas of architectural distortion estimated on a DTMslice intersecting a lesion approximately at its largest cross section;

FIG. 3 is a graph illustrating a distribution of breast density in termsof BI-RADS category for 26 breasts estimated by an MQSA radiologist fromconventional mammograms;

FIG. 4( a) is an example of a mass imaged on a DTM slice;

FIG. 4( b) is an example of a screen-film mammogram;

FIG. 5 is a flow chart illustrating a method of processing a set of DTMimages for one or more patients to screen for lesion on DTM mammograms;

FIG. 6( a) is a DTM slice intersecting a spiculated mass;

FIG. 6( b) is the segmented mass from FIG. 6( a) in the correspondingslice after 3D region growing segmentation;

FIG. 6( c) illustrates a rubber-band straightening transform (RBST)image of a 60-pixel wide region around the mass from FIG. 6( a);

FIG. 6( d) illustrates a gradient image obtained by Sobel filtering ofthe RBST image from FIG. 6( c);

FIG. 7 is an ROC curve of a linear discriminant classifier obtained fromleave-one-case-out testing;

FIG. 8 is an FROC curve for test performance of a detection system forDTM mammograms.

DETAILED DESCRIPTION

Referring to FIG. 1, a computer aided diagnosis (CAD) system 20 that maybe used to detect and diagnose breast cancers includes a computer 22having a processor 24 and a memory 26 therein and having a displayscreen 27 associated therewith. As illustrated in an expanded view ofthe memory 26, a breast cancer detection and diagnostic system 28 in theform of, for example, a program written in computer implementableinstructions or code, is stored in the memory 26 and is adapted to beexecuted on the processor 24 to perform processing on one or more setsof digital tomosynthesis mammography (DTM) images 30, which may also bestored in the computer memory 26. The DTM images 30 may include DTMimages for any number of patients and may be entered into or deliveredto the system 20 using any desired importation technique. Generallyspeaking, any number of sets of images 30 a, 30 b, 30 c, etc. (calledimage files) can be stored in the memory 26 wherein each of the imagefiles 30 a, 30 b, etc. includes numerous DTM scan images associated witha particular DTM scan of a particular patient. Thus, different ones ofthe images files 30 a, 30 b, etc. may be stored for different patientsor for the same patient at different times. As noted above, each of theimage files 30 a, 30 b, etc. includes a plurality of images thereincorresponding to the different slices of information collected by a DTMimaging system during a particular DTM scan of a patient. The actualnumber of stored scan images in any of the image files 30 a, 30 b, etc.will vary depending on the size of the patient, the scanning imagethickness, the type of DTM system used to produce the scanned images inthe image file, etc. While the image files 30 are illustrated as storedin the computer memory 26, they may be stored in any other memory and beaccessible to the computer 22 via any desired communication network,such as a dedicated or shared bus, a local area network (LAN), wide areanetwork (WAN), the internet, etc.

As also illustrated in FIG. 1, the breast cancer detection anddiagnostic system 28 includes a number of components or routines whichmay perform different steps or functionality in the process of analyzingone or more of the image files 30 to detect and/or diagnose breastcancers. As will be explained in more detail herein, the breast cancerdetection and diagnostic system 28 may include 3D gradient fieldanalysis routines 34, multi-resolution filtering routines 36, objectsegmentation routines 37, and feature extraction routines 38. To performthese routines 34-38, the breast cancer detection and diagnostic system28 may also include one or more three dimension image processing filters40, feature classification routines 42, classifiers 43, such as neuralnetwork analyzers, linear discriminant analyzers which use lineardiscriminant analysis routines to classify objects, support vectormachines, rule based analyzers, including standard or crisp rule basedanalyzers and fuzzy logic rule based analyzers, etc., all of which mayperform classification based on object features provided thereto. Ofcourse other image processing routines and devices may be includedwithin the system 28 as needed.

Still further, the CAD system 20 may include a set of files 50 thatstore information developed by the different routines 34-38 of thesystem 28. These files 50 may include temporary image files that aredeveloped from one or more of the DTM images within an image file 30 andobject files that identify or specify objects within the DIM images. Thefiles 50 may also include one or more object files specifying thelocation and boundaries of objects that may be considered as lesioncandidates, and object feature files specifying one or more features ofeach of these lesion candidates as determined by the feature classifyingroutines 42. Of course, other types of data may be stored in thedifferent files 50 for use by the system 28 to detect and diagnosebreast cancers from the DTM images of one or more of the image files 30.

Still further, the breast cancer detection and diagnostic system 28 mayinclude a display program or routine 52 that provides one or moredisplays to a user, such as a radiologist, via, for example, the screen27. Of course, the display routine 52 could provide a display of anydesired information to a user via any other output device, such as aprinter, via a personal data assistant (PDA) using wireless technology,etc.

During operation, the breast cancer detection and diagnostic system 28operates on a specified one or ones of the image files 30 a, 30 b, etc.to detect and diagnose breast cancer lesions associated with theselected image file. After performing the detection and diagnosticfunctions, which will be described in more detail below, the system 28may provide a display to a user, such as a radiologist, via the screen27 or any other output mechanism, connected to or associated with thecomputer 22 indicating the results of the breast cancer detection anddiagnostic process. Of course, the CAD system 20 may use any desiredtype of computer hardware and software, using any desired input andoutput devices to obtain DTM images and display information to a userand may take on any desired form other than that specificallyillustrated in FIG. 1.

Computer-aided diagnosis (CAD) has been shown to improve breast cancerdetection and characterization in mammography. Although a preliminaryevaluation indicated that the breast lesions can be visualized in DTMimages more easily than on conventional mammograms, the overalldetection sensitivity and specificity of DTM in comparison withconventional mammograms remain to be investigated. With DTM, the numberof reconstructed slices for each breast is very large. Even at 1-mmslice thickness, the number of slices per breast will range from about30 to over 80. The time required for interpretation of a DTM case can beexpected to be much greater than that for conventional mammograms. Withthe increase in radiologists' workload, the chance for oversight ofsubtle lesions may not be negligible. CAD will likely play a role in theinterpretation of DTM as for conventional mammograms.

The DTM images may be acquired by a GE digital tomosynthesis mammographysystem. The system may have a flat panel CsI/a:Si detector with a pixelsize of 0.1 mm×0.1 mm. The DTM system acquires projection views (PVs) ofthe compressed breast over a 50-degree arc in the mediolateral oblique(MLO) view. A total dose for the PVs was designed to be less than 1.5times that of a single standard film mammogram. DTM slices arereconstructed at 1-mm slice spacing using an iterativemaximum-likelihood algorithm.

Although a specific example of a DTM system is described here, DTMimages may be acquired by DTM systems of other manufacturers and used asinput to the CAD system 20. The digital detector can be made of othermaterials such as a selenium detector or a storage phosphor detector.The pixels of the detector can be any sizes acceptable for mammographicimaging, which may range from 0.05 mm×0.05 mm to 0.2 mm×0.2 mm. The PVsmay be acquired with different parameters such as the angular range,number of PVs per breast, radiation dose, mammographic views. The DTMslices may be reconstructed at different pixel sizes or slice thickness.DTM slices can be reconstructed using a variety of reconstructionalgorithms including, but not limited to, iterative maximum-likelihoodalgorithms, filtered backprojection algorithms, backprojectionalgorithms, algebraic reconstruction techniques (ART), simultaneous ART,and simultaneous iterative reconstruction techniques (SIRT).

A case may consist DTM slices of a single breast. In a preliminarystudy, 26 cases including 23 masses and 3 areas of architecturaldistortion were available. Thirteen of the masses and 2 of thearchitectural distortion were proven to be malignant by biopsy. The truelocation of the mass or the architectural distortion in each case wasidentified by a Mammography Quality Standards Act (MQSA) radiologistbased on the diagnostic information. The longest diameter of the lesionsranged from 5.4 mm to 29.4 mm (mean=14.2 mm, median=12.1 mm) asestimated on the DTM slice intersecting the lesion approximately at itslargest cross section. The distribution of the longest diameter of themasses or the areas of architectural distortion is shown in FIG. 2. Thedistribution of the breast density in terms of BI-RADS category for the26 breasts as estimated by an MQSA radiologist from viewing theconventional screen-film mammograms is shown in FIG. 3.

An example of a DTM slice intersecting a spiculated mass is shown inFIG. 4( a). The same mass imaged in a conventional screen-film mammogramof the same view is shown in FIG. 4( b) for comparison. The spicules ofthe mass are much more conspicuous in the DIM slice than in themammogram, likely attributable to the reduced structured background inthe DTM image.

Computerized Detection

FIG. 5 depicts a flow chart 60 that illustrates a general method ofperforming breast cancer detection and diagnosis for a patient based ona set of DTM images for the patient as well as a method of determiningwhether the detected breast lesions are benign or malignant. The flowchart 60 of FIG. 5 may generally be implemented by software or firmwareas the breast cancer detection and diagnostic system 28 of FIG. 1 if sodesired. Generally speaking, the method of detecting breast cancerdepicted by the flow chart 60 includes a series of steps 62-74 that areperformed on each of the DTM images for a particular image file 30 of apatient to identify and classify the areas of interest on the DTMimages.

The breast cancer detection CAD system 28 includes several majorsteps—prescreening, segmentation, feature extraction, false-positive(FP) reduction, as shown in the schematic in FIG. 5. For a given case,DTM slices containing the entire breast volume are input into the systemfor processing (block 62).

Prescreening:

In the prescreening step, the goal is to detect lesion candidates (block64) by using the property that the breast lesions are generally denserand have higher gradient at the transition between the breast tissuebackground and the lesion. To achieve this goal, a three-dimensional(3D) gradient field analysis can be applied to the volumetric data setof each case as described below. Of course, variations of the method orthe implementation to extract similar gradient convergence informationare possible. To reduce noise in the gradient calculation, the imagevoxels may first be averaged to obtain a smoothed volumetric data set.The gradient field analysis is calculated in a spherical region aboutthe size of a lesion, e.g., 5 to 10 mm in radius, centered at eachvoxel, c(i), of the breast volume. The gradient vector at each smoothedvoxel v(j) in the spherical region is computed and the direction of thegradient vector is projected to the radial direction from the centralvoxel c(i) to the voxel v(j). The average gradient direction over aspherical shell of voxels at a radius, R(k), of k voxels from c(i) iscalculated as the mean of the gradient directions over voxels on threeadjacent spherical shells, e.g., R(k−1), R(k), and R(k+1). Finally, thegradient field convergence at c(i) is determined by computing somestatistics, e.g., the maximum of the average gradient directions amongall shells in the spherical region. The gradient field convergencecalculation is performed over all voxels in the breast region, resultingin a 3D gradient field image.

An alternative method to detect lesion candidates is to apply 3D Hessianmatrix analysis to the DTM volumetric data set. Hessian matrix has beendescribed in the literature but has not been applied to DTM for lesiondetection. Hessian matrix is calculated by convolution of the secondderivatives of 3D Gaussian filters with the DTM data set. MultiscaleGaussian filters covering the size range of the lesions of interest willbe applied. For each scale, a response function designed to enhance 3Dspherical objects within a size range is calculated at each voxel. Theresponse functions at different scales are combined by a maximumoperation. The resulting response image depicts locations of potentiallesions as locations of relatively high responses. The response thatexceeds a selected threshold can then be labeled as potential lesions.These locations may be subjected to subsequent analysis, similar to thatapplied to the locations obtained by gradient field analysis, asdescribed below. Of course, variations in the implementation of theHessian matrix analysis are possible.

Segmentation:

Still referring to FIG. 5, a volume of interest (VOI) (e.g., 256×256×256voxels) is then identified with its center placed at each location ofhigh gradient convergence or high response. The object in each VOI issegmented by a 3D region growing method in which the location of highgradient convergence is used as the starting point and the object isallowed to grow across multiple slices (block 66). Region growing may beguided by the radial gradient magnitude. The growth of the object isterminated where the radial gradient reaches a threshold valueadaptively selected for the local object. After region growing, allconnected voxels constituting the object are labeled.

An alternative segmentation method is to use k-means clustering to groupthe voxels in the VOI into a lesion class and a background class. Thelargest connected object in the lesion class will be used as the initialregion for a 3D active contour model or level set segmentation method.The search for object boundary by the 3D active contour model is guidedby minimization of an energy function, which is a weighted sum ofinternal energy terms such as homogeneity energy, continuity energy,two-dimensional (2D) curvature energy, and 3D curvature energy, andexternal energy terms such as 2D gradient energy, 3D gradient energy,and balloon energy. The weights of the 3D active contour are optimizedusing a training set of DTM lesions.

Feature Extraction:

The 3D object characteristics are then extracted from the object (block70). Three groups of features—morphological features, gray levelfeatures, and texture features—are extracted from the segmented object.The morphological features describe the shape of the object. Theyinclude, but are not limited to, the volume in terms of the number ofvoxels in the object, volume change before and after 3D morphologicalopening by a spherical element with a radius of several voxels, e.g., 5voxels, the surface area, the maximum perimeter of the segmented objectamong all slices intersecting the object, and the longest diameter ofthe object. The compactness of the object is described in terms of thepercentage overlap with a sphere of the same volume centered at thecentroid of the object. The gray level features include, but are notlimited to, the contrast of the object relative to the surroundingbackground, the minimum and the maximum gray levels, and thecharacteristics derived from the gray level histogram of the object suchas the skewness, kurtosis, energy, and the entropy.

The texture features can be extracted either in 2D or 3D. One example of2D texture feature extraction is described as follows. On each slice,the cross section of the 3D object is treated as an object in a 2Dimage. The rubber-band straightening transform (RBST) previouslydeveloped for analysis of masses on 2D mammograms is applied to theobject. A region around the object margin, e.g., a band of 60 pixelswide, is transformed to a rectangular coordinate system. In the RBSTimage, the lesion boundary will be parallel to the long side, and theradially oriented spicules or textures are aligned and parallel to theshort side, thus facilitating texture analysis. A gradient image of thetransformed rectangular object margin is derived from Sobel filtering.Texture features can be extracted from the run length statistics (RLS)of the gradient image in both the horizontal and vertical directions.The texture features may include, but are not limited to, short runsemphasis, long runs emphasis, gray level nonuniformity, run lengthnonuniformity and run percentage. Any other type of texture features mayalso be extracted, such as spatial gray level dependence (SGLD)textures, gray level difference statistics (GLDS) textures, and the Lawstextures. Detailed description of the RBST and the RLS and SGLD texturefeatures for mammographic lesions is well known by those of ordinaryskill in the art. For a 3D object in the DTM data set, each RLS texturefeature is averaged over the corresponding feature values over slicescontaining the segmented object.

3D texture extraction methods can be considered to be an extension of a2D method. For example, the lesion can be analyzed in three groups ofplanes: (1) along the reconstructed high-resolution DTM slices that areparallel to the detector plane, (2) a set of planes that cut through thecentroid and the south-north poles, and (3) a set of planes that cutthrough the centroid and the east-west poles. On each plane, the crosssection of the lesion will be similar to a lesion in 2D. The RBST isthen applied to the band of pixels around the lesion boundary totransform it to a rectangular image. The texture in the RBST image canbe further enhanced by, e.g., Sobel filtering or multi-resolutionwavelet decomposition and reconstruction. Texture measures such as theRLS, the SGLD, the GLDS, and the Laws textures can then be extractedfrom the gray level statistics of the RBST image. The texturecharacteristics will be extracted from the planes over all orientations,described above. The corresponding texture measure extracted around thelesion can be averaged over all planes in each group to obtain anaverage texture measure over all directions around the lesion. Likewise,the RBST and texture feature extraction can be applied to the band ofvoxels inside the lesion boundary for further analysis of lesioncharacteristics.

The extracted features may incorporate object characteristics,including, but not limited to, spiculation features, boundary sharpness,and shape irregularity, that can be used by a trained classifier toestimate the likelihood that a lesion is malignant or benign. This maybe a part of the output information displayed to a user.

Feature Classification and False Positive (FP) Reduction:

If a small data set is all that is available, a leave-one-case-outresampling technique can be used for training and testing theperformance of the CAD system 20. If a large data set is available, theCAD system can be trained and tested with independent data sets. Aclassifier is trained to differentiate true lesions from false positive(FP) objects (block 72). Another classifier is trained to differentiatemalignant and benign lesions. The classifier may include, but are notlimit to, a linear discriminant analysis, a neural network, a supportvector machine, rule-based classification, or fuzzy logicclassification.

The free response receiver operating characteristic (FROC) analysis maybe used to evaluate the test performance of the CAD system. A decisionthreshold is applied to the test discriminant score of each detectedobject. For an object that has a discriminant score above the threshold,the object may be compared to the true lesion location of that case. Theobject is considered to be true positive if the centroid of the truelesion marked by the radiologist falls within the volume of the object;otherwise, it is a false positive. For each decision threshold, thedetection sensitivity and the average number of FPs per case isdetermined from the entire data set. The FROC curve is generated byvarying the decision threshold over a range of values. After the CADsystem is trained and tested to be useful, the CAD system can be appliedto new patient cases as obtained in clinical practice. An appropriatedecision threshold corresponding to a preferred sensitivity andspecificity can be chosen along the FROC curve for the clinicalapplication. The CAD system will detect the suspicious lesion locationsand display these locations to the user. The CAD system will alsoestimate and display the likelihood of malignancy of the lesion to theuser as an option.

FIGS. 6( a) and 6(b) show an example of a slice through a mass in a VOIand the mass boundary obtained by 3D region growing segmentation,respectively. An example of the RBST applied to the slice of the massand the gradient image derived from the RBST image are shown in FIGS. 6(c) and 6(d). The spicules radiating from the mass are approximately inthe vertical direction and the segmented boundary of the mass istransformed to a straight line, forming the upper edge of therectangular RBST image.

For the design of the classifier for FP reduction, a stepwise or otherfeature selection procedure may select the most effective subset offeatures from the available feature pool, thus reducing thedimensionality of the feature space for the classifier. As an example,seven features may be selected from the available feature pool. The mostoften selected features are likely to include the contrast, minimum graylevel, volume change before and after 3D morphological opening, maximumperimeter, compactness, and two RLS texture features—horizontal shortruns emphasis and gray level nonuniformity. The performance of theclassifier for differentiation of true and false lesions in the featureclassification step of the CAD system 20 can be evaluated by ROCanalysis. As an example, an ROC curve obtained from a lineardiscriminant analysis (LDA) classifier with stepwise feature selectionis shown in FIG. 7. The LDA classifier was designed with the trainingsubset in each of the leave-one-case-out cycle. The trained classifierwas applied to the lesion candidates in the left-out case such that eachobject was assigned a discriminant score. The area under the ROC curvereached 0.87 with a standard deviation of 0.02.

During prescreening, most, if not all, of the masses and architecturaldistortions are typically detected. The overall test performance of thedetection system after FP reduction for a data set of 26 cases is shownas the FROC curve in FIG. 8. The system has proven to achieve asensitivity of 85% at 2.2 FPs/case and 80% at 2.0 FPs/case. Thisperformance can further be improved if it is trained with a largertraining set.

In DTM mammography, the structured background such as the densefibroglandular tissue is suppressed in the reconstructed DTM slices.However, DTM is different from computed tomography in that theoverlapping tissues are reduced but not totally eliminated. Thetomosynthesis reconstruction leaves residues of the overlapping tissueon the DTM slices. Similarly, the shadow of a lesion can be seen in mostslices even though the actual size of the lesion may only be a fractionof the breast thickness. In addition, the voxel dimension along thez-direction (i.e., the direction perpendicular to the slices) in thereconstructed slices may be several times larger than that in the x-yplane (the planes of the slices).

The boundary of an object in the z-direction is therefore not aswell-defined as that on the x-y planes. The features extracted in 3D mayhave a strong directional dependence. For example, extracted texturefeatures may be taken along the x-y planes and a 3D texture featureobtained by averaging the corresponding 2D texture values over slicescontaining the object. Alternatively, for 3D texture analysis, thetexture features may be calculated in the shell of voxels surroundingthe object or on the planes slicing through the object centroid fromdifferent directions.

Computerized Characterization

Computerized characterization of lesions on DTM images is intended as anaid to radiologists in making biopsy or follow-up recommendations inpatient management. Classification between malignant and benign lesionsmay be performed as a separate step if the radiologist prefers to havean estimate of the likelihood of malignancy (LM) of a selected lesionduring screening or diagnostic work-up. Alternatively, the computer mayprovide a malignancy rating for each suspicious lesion found by thecomputerized detection algorithm. Computerized lesion characterizationhas the general steps—preprocessing, segmentation, feature extraction,and feature classification. The specific computer vision techniquesdeveloped for characterization of breast lesions on DTMs are describedbelow. Some of the feature extraction methods and features are similarto those developed for the lesion detection system. However, the featureclassifier is designed to differentiate true lesions and false positives(normal breast tissue) in the detection system, whereas the featureclassifier is trained to differentiate malignant and benign lesions inthe characterization system.

Preprocessing:

A rectangular-prism shaped volume of interest (VOI) enclosing the lesionthat has been either been identified manually by a radiologist orautomatically detected by the computer, is extracted. The slicethickness may first be linearly interpolated to the same as the pixelsize so that the voxels are isotropic cubes. Although DTM has reducedbackground compared with regular mammograms, dense tissue adjacent to alesion is not uncommon. The background normal tissue structure isreduced by estimating the low frequency background gray level image froma shell of voxels at the periphery of the VOI. The gray level at a givenvoxel of the background image is obtained as a weighted average of thesix gray level values of the voxels at the intersection between thenormal from the voxel of interest to each of the six surfaces of theVOI. The weighting factor of each voxel value can be taken as theinverse of the length of the respective normal. The estimated backgroundimage is then subtracted from the lesion VOI.

Segmentation:

The segmentation method outlines the main body of the lesion. The lesionwithin the background-corrected VOI is segmented using eithergeneralized k-means clustering followed by a 3D active contour model ora level set segmentation method, or 3D region growing with adaptivethresholding. In the first method, a set of pixels in 3D is grouped intoa lesion class and a background class using the k-means technique.Clustering is performed iteratively until the cluster centers of theclasses stabilize. The voxels grouped into the non-background class mayor may not be connected. A 26-connectivity criterion is used todetermine the various connected objects in the 3D space and the largestconnected object in the lesion class is used as the initial region for a3D active contour model, where the weights of the energy terms in theactive contour model are optimized using a training set for the specificapplication of lesion segmentation on DTMs. In the second method, anadaptive threshold is determined as the maximum average gradient overall radial directions around the object being segmented. The thresholdcan also be made variable over the different radial directions byestimating a local threshold over a limited-angle conical section arounda given radial direction.

Feature Extraction:

Morphological, gray-level, texture, and spiculation features areextracted. 3D morphological characteristics extracted from the segmentedlesion boundary include, but are not limited to, the volume, the longestdiameter, the surface-area-to-volume ratio, the sphericity, theeccentricity of a fitted ellipsoid, the normalized radial length fromthe object centroid, and the variance of the normalized radial length.

The gray level features for lesion characterization include, but are notlimited to, the contrast of the object relative to the surroundingbackground, the minimum and the maximum gray levels, and thecharacteristics derived from the gray level histogram of the object suchas the skewness, kurtosis, energy, and the entropy.

3D texture features may be extracted from the lesion analyzed in anumber of planes, including: (1) along the reconstructed high-resolutionDTM slices that are parallel to the detector plane, (2) a set of planesthat cut through the centroid and the south-north poles, and (3) a setof planes that cut through the centroid and the east-west poles. On eachplane, the cross section of the lesion is similar to a lesion in 2D. TheRBST previously developed for analysis of masses on 2D mammograms isapplied to the object. A region around the object margin, e.g., a bandof 60 pixels wide, is transformed to a rectangular coordinate system. Inthe RBST image, the lesion boundary is parallel to the long side, andthe radially oriented spicules or textures are aligned and parallel tothe short side, thus facilitating texture analysis. The texture in theRBST image can be further enhanced by, e.g., Sobel filtering ormulti-resolution wavelet decomposition and reconstruction. Texturemeasures such as the RLS, the SGLD matrices, the GLDS, and the Lawstextures are extracted from gray level statistics of the RBST image. Thetexture characteristics may be extracted from the planes over allorientations, described above. The corresponding texture measureextracted around the lesion can be averaged over all planes in eachgroup to obtain an average texture measure over all directions aroundthe lesion.

An alternative to extract 3D texture features is to generalize the RBSTto 3D, i.e., transforming a shell of voxels surrounding the lesion to aslab of rectangular voxels with the surface of the lesion transformed toa flat plane. The texture measures described above can then begeneralized to 3D and extracted from the slab of voxels in all threedimensions.

Spiculation features are extracted both from the 2D planes describedabove, and also in 3D. For 2D spiculation analysis, a spiculationlikelihood map around the lesion margin is calculated by the statisticsof the gradient magnitudes and directions at each pixel. For a pixel(i_(c),j_(c)) on the mass boundary, a search region is defined as theset of all image pixels that i) lie outside the mass; ii) have apositive contrast; iii) are at a certain distance (e.g., less than 4 mm)from (i_(c),j_(c)); and iv) are within ±π/4 of the normal to the masscontour at (i_(c),j_(c)). At each image pixel (i,j) in the searchregion, the obtuse angle θ between two lines is computed, where thefirst line is defined by the gradient direction at (i,j), and the secondline joins the pixel (i,j) to the mass boundary pixel (i_(c),j_(c)). Amethod based on convolution with Gaussian derivatives may be used forcomputing the gradients. The spiculation measure at a mass boundarypixel (i_(c),j_(c)) is defined as the average value of θ in the searchregion. If the pixel (i_(c),j_(c)) lies on the path of a spiculation,then θ will be close to π/2 whenever the image pixel (i,j) is on thespiculation, and hence the mean of the spiculation measure will be high.The spiculation measure may be computed for a number of contours, whereeach contour is obtained by expanding the previous contour by one pixelat all pixels on the contour, and the first contour is given by thesegmentation method. The resulting image in a band around the mass isreferred to as the spiculation likelihood map. High values in the mapindicate a high likelihood of the presence of spiculations. A thresholdabove which the voxel is considered to be part of a spicule isdetermined by training with case samples. The features in this map areused to classify lesions as spiculated and non-spiculated. Spiculationmeasures extracted from the spiculation likelihood image include, butare not limited to, the number, the density, the total area ofspiculations, and the strength of the spiculation gradients to measurethe degree of spiculation of the lesion. After extracting thespiculation measures from each plane, the corresponding measures can beaveraged over all planes in the same group to derive 3D spiculationmeasures. These spiculation measures can be used directly, and incombination with morphological and texture features, as input featuresto the malignant-benign classifiers. Additional features fordistinguishing fibrous tissue overlapping with the lesion from truespiculations include measures for gradient strengths and directioninside the lesion, and the estimate of the direction and the continuityof the structures from the inside to the outside of the lesion.

For 3D spiculation analysis, the spiculation likelihood map isgeneralized to 3D by calculating the statistics of the gradientmagnitudes and directions in 3D space at each voxel in the lesion marginregion. The gradient statistics at a point are calculated over acone-shape region around the normal at that point. High values ofgradient statistics in the lesion margin indicate possible presence ofspiculations. A threshold above which the voxel is considered to be partof a spicule is determined by training with case samples. 3D spiculationmeasures such as number, density, and the proportion of voxelsconsidered to be spicules relative to the total number of voxels in thelesion margin are derived from the 3D spiculation likelihood map. Thesefeatures can then be used directly, or in combination with othermorphological and texture features, as input predictor variables ofclassifiers for differentiation of malignant and benign lesions.

Feature Classification:

A classifier is trained to differentiate malignant and benign lesionsbased on the extracted features above. The classifier may include alinear discriminant analysis, a neural network, a support vectormachine, rule-based classification, or fuzzy logic classification. Forthe design of the classifier, a stepwise or other feature selectionprocedure may select the most effective subset of features from theavailable feature pool, thus reducing the dimensionality of the featurespace for the classifier.

The formulation of the classifiers is well known by those of ordinaryskill in the art, but the combinations of segmentation and featureextraction techniques that are developed for the DTM images make thedesigned classifier specific for the CAD task described herein.

For DTM imaging, the raw data are acquired as projection view (PV)mammograms. The number of PVs will depend on the design of the DTM x-raysystem. For example, a first generation GE prototype DTM system acquired11 PVs. A second generation GE prototype DTM system acquired 21 PVs. Thetotal dose of all PVs used about 1.5 times of the dose of a conventionalmammogram. A PV is therefore noisier than a conventional mammogram.However, the PVs offer the advantage that a lesion will be projected ata slightly different angle and thus having somewhat differentoverlapping tissues on each view. A lesion that may be camouflaged bydense tissues on some views may become more conspicuous on other views.In addition, overlapping tissues that mimic lesions on some views may beless lesion-mimicking on other views. If a CAD system for lesiondetection is applied to the PVs, the complementary information on thedifferent PVs may be utilized to improve sensitivity and reduce FPs.

Furthermore, although the disclosed algorithm uses DTM slicesreconstructed from 11 PVs using the iterative maximum-likelihoodalgorithm as input, image processing methods should not depend stronglyon the reconstruction method or the number of PVs for generating the DTMslices as long as the image quality of the reconstructed slices arereasonable. The raw PV mammograms may also be directly used as input tothe CAD system, with a chosen reconstruction method installed as apreprocessing step in the CAD system 20. Alternatively, the CAD systemmay first process the PV mammograms individually using 2D imageprocessing techniques. The detected objects on the PVs may then bemerged based on the extracted features. The merged information is usedto differentiate true lesions or false positives.

When implemented, any of the software described herein may be stored inany computer readable memory such as on a magnetic disk, an opticaldisk, or other storage medium, in a RAM or ROM of a computer orprocessor, etc. Likewise, this software may be delivered to a user or acomputer using any known or desired delivery method including, forexample, on a computer readable disk or other transportable computerstorage mechanism or over a communication channel such as a telephoneline, the Internet, the World Wide Web, any other local area network,wide area network, or wireless network, etc. (which delivery is viewedas being the same as or interchangeable with providing such software viaa transportable storage medium). Furthermore, this software may beprovided directly without modulation or encryption or may be modulatedand/or encrypted using any suitable modulation carrier wave and/orencryption technique before being transmitted over a communicationchannel.

While the present invention has been described with reference tospecific examples, which are intended to be illustrative only and not tobe limiting of the invention, it will be apparent to those of ordinaryskill in the art that changes, additions or deletions may be made to thedisclosed embodiments without departing from the spirit and scope of theinvention

1. A method for using computer-aided diagnosis (CAD) for digitaltomosynthesis mammograms (DTM), comprising: retrieving a DTM image filehaving a plurality of DTM image slices; applying a three-dimensionalanalysis to the DTM image file to detect lesion candidates; identifyinga volume of interest and locating its center; segmenting the volume ofinterest by a three dimensional method; extracting one or more objectcharacteristics from the object corresponding to the volume of interest;and determining if the object corresponding to the volume of interest isa breast lesion or normal breast tissue.